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ABSTRACT 

Using A^-body simulations {N ~ 10^ — 10^), we examine how a non-axisymmetric dark 
halo affects the dynamical evolution of the structure in collisionless (stellar) discs. We 
demonstrate how the model parameters such as mass of the halo, initial conditions 
in the disc and the halo axes ratio affect morphology and kinematics of the stellar 
discs. We show that a non-axisymmetric halo can generate a large-scale spiral density 
pattern in the embedded stellar disc. The pattern is observed in the disc for many 
periods of its revolution, even if the disc is gravitationally over-stable. The growth 
of the spiral arms is not accompanied by significant dynamical heating of the disc, 
irrelevant to its initial parameters. We also investigate transformation of the dark 
halo's shape driven by the long-lived spiral pattern in the disc . We show that the 
analysis of the velocity field in the stellar disc and in the spiral pattern gives us a 
possibility to figure out the spatial orientation of the triaxial-shaped dark halo and to 
measure the triaxiality. 

Key words: Galaxies: evolution; galaxies: haloes; galaxies: kinematics and dynamics; 
galaxies: spiral 



1 INTRODUCTION 

Both dynamic models of real galaxies and cosmological sim- 
ulations give a strong evidence that the stellar discs of spiral 
galaxies are embedded into massive dark haloes. The halo 
mass and its influence on the disc evolution is the matter of 
debates. Among the different ways to put constraints on the 
halo mass it is worth mention the use of gravitational sta- 
bility criteria in the disc's plane (ISafronov I196QI: iToomre 



Lvnden-Bell and Ostrikei] Il967l : I Khoper skov et al 
Zasov et all l2QQ4l \ a condition of stellar disc 



stability 

ISotnikova et al.1 [20a 



to bending perturbations (IZasov et al.l Il99ll : 
Khoperskov et al.1 120101)" simulation 



of polar rings (Whitmore et al.l I l987: Sac kett and Sparkj 
I1990I : iReshetniko^ and SotnikQ^r2000: lodice et al.ll2Q03l) 
and the gravitational lensing ([Koopmans et all Il998l : 
iMaller et alJl2000l : lOguri et al.ll2003f L 

Morphology of galactic spiral patterns also indicates 
significant amount of the dark matter within the opti- 
cal radius jA thanaso ula et al.l Il987l : iKorchagin et al.l |200QI : 
iGriv and Ge dalin 20^). Estimates of the dark mass con- 
tribution can also be obtaine d from t h e vert i cal structure 
of ste ll ar discs analysis (e.g. iBahcall' ('l984'); ^ Zaso v et al. 
1991 ): iBizvaev and Mitronova (2002); Zasov and Bizvaev 
2002); IKorchagin et al 1 ([20031 ): IBizvaev and Mitronoval 



(|2009l ): IJust and Jahreil (|2010l ): lGarbari et al.1 (|201ll )). Max- 
imum disc models suggest that Mh/{Md + Mh) ^ 0.5 
1 . where Mh is the halo mass within the optical ra- 
dius, Mh is the mass of stellar bulge, Ma is the mass 
of stellar disc (fBottemal [1993). Note that some galac- 
tic discs (especially low luminous ones) are less mas- 
sive than that required by the maximum dis c assump- 
tion jZasov et al.1 (I20q3): iBershadv et al.1 (|201ll ): IWestfalll 
(I2OIII ): ISaburova and Zasovl (l2012l )T At any case, one may 
conclude that in most galaxies (at least of normal surface 
brightness) masses Mh and Md are comparable within the 
optical borders. 

Modern cosmological experiments sugges t a triaxial dis- 
tribut io n of the dark matter i n the haloes (e.g.lZentner et al.l 
(I2OO5I): iKuhlen et al.1 (|2007l ). iMuoz-Cuartas et al.1 (|201ll ). 
iKlvpin et al .1 (I2OIII ). and references therein). Intriguing 
problem is a possibility to explain the generation and sup- 
port for the long-lived and extended spiral structure in the 
stellar discs by triaxial potential of the dark halo. True mech- 
anism of the formation of the spiral structure is still a subject 
of discussions. Although a bar or a close satellite may induce 
the spiral density wave, these events provide a relatively 
short term support. In contrast, the dark halo can keep its 
triaxiality and support the spiral structure over much longer 
time scale. The problem is to estimate whether the existing 
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Figure 1. Evolution of the stellar velocity dispersion Cr at dif- 
ferent radii r in a model with spherically symmetric halo. The 
curve 1 demonstrates the central zone of the disc, the curve 8 - 
the periphery of the disc. The arrow indicates the time interval 
corresponding to one period of rotation tq at the disc's periphery. 



deviation from the axial symmetry is quantitatively enough 
to support the spiral pattern in collisionless gravitationally 
stable discs. 

Usually the distribution of the volume density in the 
halo is approximated as follows: 



Qh{x,y,z) = QhiO , ^ 



g2 ^2 ^2 ' 



(1) 



where the a — h — c case corresponds to the spherical sym- 
metry in the density distribution. The model with a = 6 / c 
corresponds to an axisymmetric halo coplanar with the stel- 
lar disc. Note that in this labeling convention, the relative 
ordering of a, 6, and c is unconstrained. The ratios q = b/a^ 
s — c/a describe the deviation of the density profiles front 
the spherical symmetry. For simplicity, our further consid- 
eration of the dark matter (DM) haloes assumes that a ^ b 
and c = a, where a and b are the halo major and minor axes 
in the disc plane. 

Cosmological A/'-body simulations show that the tri- 
axiality parameters s and q depend on the halo mass 
and distance to the centre dBailin and Steinmet"3 l2005l : 
iHavashi et alll2007l : Iwidrowl l2008l ) . Simulations conducted 
with a large number of particles (as the Via Lactea 
with N = (1 - 2.3) • 10^. Ipiemand. Kuhlen and Madaul 
(I2OO7I): [Kuhlen et al.l (l2007l ): IPebattista et al.1 (|2008l )T 
iDiemand and Moord (120091 )7 reveal that the dark haloes are 
more asymmetric at their central regions than at the periph- 
ery where the ratios b/a and c/a are close to to unit. For 
the topic of this paper the most important parameter of the 
halo's shape is the axes ratio q. 

Observing constraints on the dark halo geometrical pa- 
rameters come from studies of globular clusters, variable 
stars (such as RR Lyrae), blue horizontal branch stars, low- 
met allicity_gmiits^_^ed giants, Cepheids (W Vir type) and 
F-dwarfs ([Keller et al.ll2008l V Anisotropy of the velocity dis- 
persion of objects in the stellar halo of our Galaxy indicate s 
an asymmetric distribution of DM ([Morrison et al.ll2009l V 
Stellar haloes in other spiral galaxies also demonstrate sig- 
nificant deviation from sp herical symmetry in the haloes: 
c/a ^ 0.3 is suggested by ICooper et all tOld ). Accordine: 
to the Sloan Digital Sky Survey (SDSS) and Sloan Exten- 
sion for Galactic Understanding and Exploration results. 



our galac tic spheroid has a t riaxial shape with s = 0.63, 
q = 0.75 (iNewberg et al.ll2007l ). 

The polar ring galaxies are un ique obje cts for investiga- 
tion o f the DM halo shape ([Combe s 2006). Iwhitmore et al.l 
([19871 ) analysed the polar ring kinematics in several galaxies 
and found significant deviation of the DM distribution from 
the spherical sha pe: c/a ^ 0.5 (ISackett and Pogg3 1 19951 ) , 
c/a - 0.3 - 0.4 (Sack ett et allll994l ). A prolate shape o f 
the dark matter haloes was inferred by llodice et all (l2003l ) 
based on the analysis of location of the polar ring galaxies 
on the Tully-Fisher diagram. 

Dynamics of decaying satellites also helps to constrain 
the dark hal o properties, such as its shape, orientation, mass 



Helmill2004l:lFellhauer et al.1 



distribution (llbata et al]l2001bl: 

I2OO6I ). According to ISesar et al.l ([201 ll ). SDSS data suggest 
that the Milky Way's dark halo is oblate (s = 0.7) and non- 
axisymmetric in the disc plane [(q = 0.98) within 28 kpc of 
the centre. 

An assumption of the axissymmetric potential for the 
Milky Way halo applied to the Sgr dSph tidal stream yields 
controversial results: while the leading arm of the stream re- 
quires an oblate halo, the radial velocity distribu tion of stars 
in th e arm agrees better with a prolate halo ([ Vivas et al.l 
I2OO5I ) . Note that the intro duction of a triaxial potential helps 
to resolve this ambiguity ([Law et al ][2009). 

Statistical analysis of shapes of the spiral galaxies 
also reveals some deviation from the disc axial symmetry, 
whi ch in turn ma y evidence the halo triaxiality. According 
to iRvdenI (I2OO6I ). the disc ellipticity is 0.02 - 0.3 in the 
Kg photometric band. The ellipticity varies in dependence 
of th e photometri c band and of the galaxy morphological 
type ([Rvdenll2006l ). 

Despite direct and non-direct evidences of the halo 
triaxiality mentioned above, its shape and ellipticity in 
specific galaxies remains poorly known. Understanding 
of the role that the asymmetry of the halo poten- 
tial can play in dynamical e volution of stellar d iscs i n 
galaxies remains i ncomp lete. El-Zant and Hableil ( 1998 ): 
Ideta a nd Hozumi (2009); 'Tutukov and Fedoroval (120061): 

Berentzen and Shlosra an (2006); Havashi an d Navarrol 

(I2OO6I ): iHavashi et al.1 | 20071 ) studied different aspects of 
how a non- axisymmetric gravitational potential affects 
the galactic disc's dynamics. The triaxial halo allows to 
explain the observed shape of rotati on curves in low s urface 
brightness (L SB) galaxies (|Havashi and Navarrol I2OO6I : 
IWidrowlliooi ). 

Galaxies with warped discs were successfully repro- 
duced withln__th^_ji^^ 

halo ([Dubinski and Chakrabartvll2009l : iRoskar et al]l2010l ). 
The effects of the bar response to the triaxial halo perturba- 
tions and bar dissolution in prola te hal o were studied by 
iMac hado and Athanassoulal ([2OIOI ) and iKazantzidis et all 
( 2016 ). Numerous A^-body simulations are usually used 
to study the formation and evolution of stellar bars em- 
bedded in triax i al da r k matter haloes (I Athanasoulal I2OO2I : 
iBerentzen et al.1 I2OO6I : iDubinski et al.1 l2009l ). The triaxial 
halo can destroy stellar bars due to the angular momentum 
exchange (Berentzen and Shlosman 2006). The decreasing 
of the bar pattern speed was demonstrated in model s with 
tria xial dark matter halo ([Ceverino and Klvpinll2007l ). 

iBekki and FreemanI ([20021 ) explained the spiral struc- 
ture extended far beyond optical limits in a blue compact 
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Figure 2. The surface density isolines a(r, (f) at different time 
t in the model with £ = 0.1 (cr is shown in the log scale). One 
period of rotation at the disc's periphery corresponds to about 3 
our dimensionless units. 



dwar f g:alaxy NGC 2915 b y the presence of a triaxial dark 
halo. iMasset and Bureaul (|2003h modeled the spiral struc- 
ture in the outer region of the gaseous disc of NGC 2915 
generated by a triaxial dark matter potential. The main re- 
sult of these papers is that triaxial DM halo is able to pro- 
duce a two-ar ni spiral pattern in externa l regions beyond the 
optical radius. [Tinker and RvdenI (|2002h investigated effects 
of rotating triaxial haloes in the disc galaxies using A/'-body 
simulations with the rigid halo potential and found that the 
average pitch angle of the spiral structure strongly anti- 
correlates with the total mass of the system (disc+halo). 
iDubinski and Chakrabartvl (|2009h applied A^-body simula- 
tions to demonstrate the formation of warps and spirals in 
isolated discs embedded in triaxial haloes. They showed that 
the spiral pattern may emerge in gravitational stable stel- 
lar disc through interactions with a non-axisymmetric halo. 
Note that the evolution of the grand design spiral structure 
in stellar discs embedded in the triaxial dark matter haloes 

remains poorly studied. 

In our previous paper (iKhoperskov et aLllioTj ) we stud- 
ied the dynamics of a gaseous galactic disc in a non- 
axisymmetric gravitational potential in general. In the 
present paper we perform A/'-body modelling (with N ^ 
10^ — 10^) of stellar collisionless discs embedded in a non- 
axisymmetric dark halo with focusing on the formation of 
the spiral pattern, its morphology and kinematics. In Sec- 
tion 2 we describe our modelling technique and the initial 
conditions in the our galactic discs. The spiral structure for- 
mation and the kinematics of the density waves are analysed 
in Section 3. Discussion and main conclusions are given in 
Sections 4 and 5, respectively. 



2 MODELLING COLLISIONLESS STELLAR 
DISCS 

Throughout this paper we study the dynamics of the galactic 
discs in the frames of a self-consistent A/'-body problem by 
solving the motion equations for equal mass test particles: 



E 



j = 1,2,..., AT, (2) 



here Vj = {x^y, z} is the coordinate vector of j -th particle, 
Md is the total stellar disc mass, ^^(x,^/,^;) is the bulge 
potential, t) is the halo potential and m = Md/N 

is the mass of a test particle. 

The force of gravity in the disc is determined 
by the summation over all particles in Equation (2) 
and c alculations via the TreeCode Top Down algo- 
rithm ([Barnes and "m^ll986l ). Initial distribution functions 
in the equilibrium st ellar disc are described in de tail in 
our previous papers (iKhoperskov et al.1 l200l l2010l ). The 
schemes of the 2nd and 4th orders we apply provide a good 
accuracy of integration of Equations (2). 

We introduce the following dimensionless model units: 
Our spatial distance unit is assumed to be 10 kpc, the grav- 
itational constant is chosen to be G = 1, the time unit is 91 
Myr, the velocity is expressed in the units of 110 km s~^. 
Assumptions G — 1^ V — 1 lead to the total mass of model 
galaxy equal to about 10^^ M© within R= 1. 

The number of particles in our models is A = (1 — 10) • 
10^. The potential for i-th particle at j-th point is ^fij — 

—vfij ^ rf- + , where the cut off radius is Vc ^ 10~^ (which 

corresponds to ^ 10 pc). The opening angle parameter G is 
assumed to be 0.1 — 0.5, although we did not find significant 
variations of the results between models with different G. 
The results presented in this paper were obtained assuming 
G = 0.1, if alternative value is not specified. This values 
of Tc and G determine the accuracy of the calculations. The 
integration step is chosen to be 2-10~^, which approximately 
corresponds to 2 • 10^ yr. 

Basic model parameters are constrained in the range 
suggested by observations as follows: 

• Relative dark halo mass 1 ^ fi = Mh/Md ^ 4 within 
the limits of the stellar disc R = Avd = 1, where is the 
radial scalelength of the disc. 

• lih — Mh/Md = — 0.3 is the relative stellar mass of 
the bulge with scalelength = 0.1 — 0.3. 

• To ~ 3 is the period of rotation at the periphery of the 
disc. 

• Nonaxissymmetric halo parameter e — 1—b/a = 0—0.25 
in the disc plane (3- 

• = — 0.3 is the angular rotation velocity of the halo. 
The halo corotation radius is located far outside the stellar 
disc (Qh = 0.3 or about 3 km s~^ kpc~^). 

• Th, = — 20 is the time scale of variation of the non- 
axisymmetry parameter from to £. 



^ It should be emphasised that the value £ = 1 — b/a calcu- 
lated for the dark halo volume density distribution differs from 
that for the gravitational po tential e^. Approximately, ~ Se 
(|Binnev and Tremain"3l2008l V 
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Figure 3. Evolution of the disc parameters averaged over the azimuthal angle at various distances from the centre (curve number 
increases from the centre towards the periphery). Strong fluctuations at the beginning of the experiment {t < 5) are due to the deviation 
from precise balance at the initial moment. 



We assume the exponential surface density distribution 
in the disc as an initial condition: a{r) = ctq exp(— r/r^), 
where = 0.25 of dimensionless distance units is the radial 
scale length in the disc. At the initial stage, all stellar disc 
particles are situated within the radius r ^ 5rd — 1.25. 
We assume a flat rotation curve with the maximum velocity 
V^max — 2 in dimensionless units, which corresponds to the 
rotation period of ^ 3 at th e disc's periphery. 

Followine; iBegeman" et al.1 (|l99lh . we use a pseudo- 
isothermal model of the density distribution in the halo: 

QhQ 



Qh{r) = 



(3) 



l + (r/a)2 • 

This distribution of the dark mass provides a constant cir- 
cular velocity Vc for r ^ a. The gravitational potential can 
be written 



ar ctg(0 



^1. 1 + r 



(4) 



([Khoperskov et al]l2012l ). which corresponds to the density 
law ([3)) for the case of a = b = c. 

Note that the cosmological simulations predict that viri- 
alised dark matter haloes should have a universal density 
profile. For a spherically average d density distribution it can 
be written (jNavarro et al.lll997l l: 
Qhi 



Qh{r) 



(5) 



{r/a){l + {r/a)r 

Its distinctive feature is the inner density cusp, which 
however may be dissolved due to its interaction with 
baryonic matter in the central region of g alaxie s 
(see f.e. ICu stafsso n et al.1 (l20Q6l ): iMaccio et all (l2012l ): 
I Khoperskov e t al. (2 012l lV Moreover, inevitable singularity 
emerges at r = in the case of transformation from a 
spherical NFW profile to the triaxial shape. Another widely 
used dark matter density distribution is Burkert's pro- 
file (IBu rkert 1995), which usually fits observing data very 
well ([Centile et al. 200^). It is very close to the NFW profile: 
it has a cored density distribution in the inner part (r < a). 
Several more halo dens i ty pro files were introduced, such as 
expo nential profile IFuxI (Il997l ) and Einasto profile (|Einastol 
Il965l ). All density distribution functions mentioned above 
(except ([5])) allow generalization of the spherical dark mat- 
ter distribution shape to the triaxial case using the formal 
replacement r ^ ^ (see eq.[T]). Note however, that the choice 
of the accepted law for Qhir) makes only a small effect on 
the results described below. 



Rotation pattern in the dark haloes predicted by 
AC DM cosmological models is similar to that in giant el- 
liptical galaxies, where the kinetic energy of the rotation is 
much smaller than that in the stellar chaotic motions. The 
origin and properties of angular momentum of DM haloes 
were in vestigated by IVitvitska et aP (|2002l ): iMaccio et al.l 
(l2QQ7l ): I Bett et a ] (l2007l ). For the disc galaxies cosmologi- 
cal simulations suggest slow rotation of their dark matter 
haloes within Ropt- The rotation (tumbling) of the over- 
all shape of the halo in ACDM simulations was d is cusse d 
by iBailin and Steinmetd (|2QQ4l ): iBrvan and Cress I (|20Q7l ). 
Authors conclude that the tumbling rates of the haloes are 
less than one km s~^ kpc~^. Evidently, a disc subsystem 
rotates much faster than a halo. Hereafter we use a non- 
rotating or slow rigidly rotating model halo with the angular 
velocity Qh, which places the corotation radius at the disc's 
periphery Qh ^ Vma^/ Ropt- Our principal results are shown 
for the case of rigid DM halo, while the situation with the 
live halo is considered in Section 3.4. 



3 DISC DYNAMICS: PROPERTIES OF THE 
SPIRAL STRUCTURE 

Propagation of the density waves of gravitational nature 
in a self- gravitating stellar disc requires that the stel- 
lar velocity dispersion were close to the stability condi- 
tio n Qt ^ 2 (or Qt ^ 3 at the periphery of the disc. 



see iKhoperskov et al.l (|20Q3 |)). where Qt is the Toomre sta- 
bility parameter ( Toomrd 1 196 4). The case Qt = 1 is un- 
stable due to non-radial perturbations. In the following sec- 
tion we consider gravitationally stable discs (Qt ^ 2 — 4 at 
r (1 — 2) • r^) and study the pure effects from the non- 
axisymmetric dark matter distribution. We also compare the 
properties of the spiral pattern generated by triaxial haloes 
in the cases of stable and unstable discs (see Section 3.3). 



3.1 Dynamically stable stellar discs 

The initial distributions of the radial Cr, azimuthal 
and vertical velocities dispersions Cz are chosen to be high 
enough to prevent the disc gravitational instability, which is 
able to generate spiral waves in the initially axisymmetric 
model by its own. The stellar velocity dispersion should be 
high enough to support the gravitational stability against 
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Figure 4. Amplitudes of the Fourier harmonics for m = 2 mode 
for the entire disc in different models: 1 : /j. = 4, 2 : /j. = 3, 3: 
II = 2, 4'- = ^ = 0.15, rij = 0.3 for all models). 
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Figure 5. Mean amplitude of the time- averaged Fourier harmon- 
ics for the entire disc in the two-arm mode as a function of the 
parameter s. 

perturbat ions in the disc plane an d against the bending in- 
stability (|Khoperskov et aP |201Q| ) . An iterative procedure 
is applied fo r making marginally stab le 3-D stellar discs as 
described bv lKhoperskov et alJ (|2003l ). 

First, we tested the disc's stability in the potential of the 
spherical halo. We assureed that the initial disc parameters 
(cr(r), Cc^(r), Cz{r), rotation velocity V{r), surface density 
cr(r), and the disc vertical scale h{r)) do not change signifi- 
cantly during several revolutions of the disc. As an example, 
figure 1 illustrates the temporal variation of the Cr(r). Col- 
lisionless stellar disc in this model does not undergo any 
dynamic heating neither to binary particle interactions nor 
due to collective processes. Under the chosen conditions the 
disc remains axisymmetric during dozens revolutions (which 
corresponds to a few billion years). 

3.2 Generation of the spiral structure in the 
stellar disc 

Our simulations started from the spherical halo and axisym- 
metric stellar disc, then we slowly and adiabat ideally (over 
time scale rn 3) increased the halo's non-axisymmetric 
parameter s. 

We performed about 30 modelling runs with different 
sets of initial model parameters mentioned above. It is ev- 



ident that in all models with e > 0.05 the presence of the 
non-axisymmetric dark halo excites the formation of a two- 
arm spiral pattern in the stellar disc. In figure 2 we show 
snapshots of the typical evolution of the surface density 
in the model with e = 0.2. Within 1-2 periods of rota- 
tion To, a two- armed spiral structure forms in the initially 
axisymmetric exponential disc. A nested B-like structures 
and a central bar can be identified in the inner disc region 
(r ^ 0.2 ^ = 0.25). The rigid triaxial halo models (mod- 
els without the disc feedback to the halo) generate a long 
lasting global spiral pattern that extends from the centre to 
the disc periphery. It is important to note that interactions 
between the disc and non-axisymmetric dark halo's poten- 
tial do not introduce any significant systematic variations to 
the azimuthally-averaged disc parameters. 

Evolution of Cr(r), V{r) and a{r) is shown in figure 3. 
Note that the centre of the disc (curves 1 and 2) remains un- 
perturbed over the course of the simulations, whereas quasi- 
periodic radial oscillations can be observed, especially in the 
outer part of the disc. In contrast to the global gravitational 
stochastic s piral structure that dev elops in the slightly un- 
stable discs (|De Simone et al.ll2004h , here the heating of the 
stellar disc is insignificant (except for the short initial stage 
of the experiment). 

The model parameters /i, /i^, rb and th significantly 
affect the growth increment of the spiral structure. The am- 
plitude of the stellar density wave at early stages in the sim- 
ulations is primarily determined by the relative halo mass 
parameter /i (see figure 4). However, after several tq, two- 
arm spirals reach the similar amplitudes at a fixed value of 
£, independently of /i. Our modelling demonstrates the de- 
pendence of the eventual amplitude A2 of the spiral waves 
on the triaxiality parameter e for the spiral mode m = 2 (see 
figure 5). It seems that £ is a crucial parameter, which de- 
fines the properties of the strength of the spiral structure 
as a whole. Despite the fact that the growth rate of the 
spiral waves depends on the relative mass of the halo (see 
figure 4), the final amplitude of the waves is determined by 
the non-axisymmetric parameter e (see figure 5). 

There are two distinctive properties of the spiral struc- 
ture generated by the triaxial halo. First, the spirals are 
generated in a sufficiently hot (i.e. gravitationally stable) 
disc with the Toomre parameter Qt 2 — 4, in contrast to 
the standard density waves theory. Therefore, this mecha- 
nism can produce spirals even in the over-stable stellar discs, 
e.g. in the outer regions of galaxies far beyond their opti- 
cal size limits, where the disc surface density is much lower 
than the gravitational stability threshold. Second, the for- 
mation of the spiral structure leads to a weak stellar disc 
warming-up during the initial stage of the spiral pattern 
development, and after that the velocity dispersion compo- 
nents only slightly increases even in the case of high am- 
plitude spiral waves (figure 3). These properties are signif- 
icantly different from those seen in A/'-body simulations of 
gravitationally unstable discs, in which the development of 
the gravitational instability is accompanied by significant 
growth of the velocity dispersion during one period of rota- 
tion. It should be noted that a slight disc heating for large 
amplitude waves that we observe in the our modelling may 
be a result of limited accuracy in the numerical solution of 
the motion equations (2). Indeed, we found that the disc 
heating decreases with the increasing of the number of test 
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particles. The heating up effect may also be attributed to 
the lack of perfect disc equilibrium at the beginning of the 
simulations. The final conclusion about the low disc heating 
requires modelling with higher spatial resolution. 



3.3 Simulations with different stellar disc initial 
state 

While we started simulations with the initially stable disc in 
the previous section, we should consider whether the initial 
gravitational instability of the disc can qualitatively change 
the results of section 3.2. To reveal the effect of self- gravity 
we have performed simulations for initially unstable disc by 
assuming the relative halo mass Mt/Md — 1 and Qt — 1. 
We started with the axisymmetric halo to be sure that the 
evolution of the disc and spiral structure is fully determined 
by the global gravitational instability of the entire disc, and 
that the formation of the spiral waves is accompanied by 
the disc heating up to the marginal stability condition (see 
Figure 6). After that we ran another simulation with the 
same initial conditions, but assuming the halo with e = 0.1. 
Relatively low halo mass led to slower growth of the wave 
amplitude than that shown in Figure 4 for /i = 1. Curi- 
ously, the final amplitude of the spiral waves is practically 
the same independently on the mechanism of their genera- 
tion (either non-axisymmetric halo + gravitational instabil- 
ity, or pure gravitational instability). A general character of 
the evolution of the velocity dispersion is also very similar 
in both models. It confirms that the non-axisymmetric halo 
does not contribute essentually to the eventual disc heating. 

To analyze the direct role of a triaxial halo and to sep- 
arate it from effects of the disc initial state, below we com- 
pare different numerical models of evolution of gravitation- 
ally stable stellar discs in the halo field. It was assumed that 
the disc always remains in a quasi-equilibrium state even it is 
non-axisymmetric, since the halo shape changes slowly and 
adiabatically in our models. Model (a) is a "reference" model 
for the initially round-shaped disc described in Section 3.2. 
In model (6) we introduce a very large timescale for the halo 
triaxiality growth {th = 15) for initially round disc to be sure 
that the disc has enough time to follow the slowly changing 
halo potential. In our third model (c) we use the initially el- 
liptical disc with the same spatial major axis orientation as 
for the non-axisymmetric halo. In our last model {d) we took 
the results of long-time evolution of a "reference" model of 
round disc in the non-axisymmetric halo at t = 20. After 
that the spiral structure was radially smoothed (which in- 
troduced slight ellipticity to the disc), and next simulation 
was started from this equilibrated, spiral-free initial state. 

Results of modelling are shown in Figure 7. In all four 
models the spiral structure emerges. The main difference 
between the models is observed at the central part of the 
disc. In models (a) and (6) a G-like structure forms. In the 
model (c) of initially elliptical disc a large bar was formed, 
and beyond the bar there is a tightly wound spiral structure. 
In the case {d) there is a small bar and very even spiral arms. 
Formation of a bar in the potential of triaxial DM ha l o link s 
our models with those developed bv lBerentzen et al.1 (|2006h . 
where the bar generation and bar-halo interaction are also 
described. It is essential that our models demonstrate the 
emergence of the spiral arms in the very outer regions of the 
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Figure 6. Top panel: evolution of the Fourier harmonic m = 2 
in a gravitationally unstable disc {Qt = 1) embedded in a non- 
axisymmetric halo £ = 0.1 (solid line) and in a symmetric halo 
(dashed line). Bottom panel: evolution of the velocity dispersion 
in the same two models at the radius of r = 0.5. 




Figure 7. Disc surface density at t = 8 in models with different 
initial state of the stellar disc for the non-axisymmetrical halo 
with e = 0.1: (a) the reference model of the symmetric disc; (b) 
evolution with long term quasi adiabatic growth of £ with = 20; 
(c) the model with a non-axisymmetric disc; (d) the model started 
with a radially smoothed disc made with those from the model 
(a) after t = 20 evolution. 



disc driven by the halo triaxiality, independently of the disc 
gravitational stability and its initial shape. 

To analyze the properties of the spiral structure with- 
out considering the central disc regions, we considered the 
evolution of the Fourier harmonics and velocity dispersion 
beyond the radius r = 0.5 (the solid circle in Figure 7). De- 
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Figure 8. Analysis of the spiral arms in different models beyond 
r = 0.5 — the designations are kept the same as in figure [3?3] Top 
panel: evolution of the Fourier harmonic m = 2. Bottom panel: 
evolution of the velocity dispersion. 



spite the difference in the disc initial state, the amplitude 
and its time-dependent evolution are rather similar (see Fig- 
ure 8, top panel). The final amplitude of the spiral waves is 
about 10% in all cases. The initial values of the Fourier har- 
monics m = 2 in (c) and (d) models are large because of 
the initial ellipticity. The model (6) with a long timescale 
of the shape formation of adiabatic halo shape transforma- 
tion reveals a slow quasi-steady increasing of the spiral wave 
amplitude, which follows the slow growth of the halo triax- 
iality. These numerical experiments confirm that formation 
of the spiral structure in the models is entirely due to in- 
teraction of the disc with the non-axisymmetric halo. After 
the spirals have formed, there is no significant dynamical 
heating introduced to the disc (see bottom panel in Figure 
8). Some initial growth of the velocity dispersion observed 
in the model (c) is caused by the relaxation processes in 
the elliptical-shaped disc while it transits to another (also 
elliptical) disc state due to interaction with equally oriented 
halo. After this transition has finished the disc heating is 
ceases. 



3.4 Kinematics of the spiral structure 

The non- stationary evolution of the spiral waves is charac- 
terised by quasi-periodic oscillations that are caused by the 
radial propagation of spiral waves (figure 3) . The spiral wave 
pattern rotates with the angular velocity Qp, which slightly 
changes with time. Figure 9 shows the evolution of the wave 
amplitude at different radii. Estimation of the phase shift 
along the radius enables us to find the pitch angle of the 
spiral pattern, which ranges from 5° to 30°. For illustration, 
we assume it to be ^ 12° in our model with s — 0.1. In our 
simulations, the pitch angle is a function of radius: it grow 
from the centre to the periphery in the galactic disc, so the 
spiral arms become more open. This radial dependence may 
be used for the revealing of the halo non-axisymetry signa- 
ture in the real galaxies. 




Figure 9. Evolution of the amplitude of the Fourier harmonic 
m = 2 in the surface density distribution at various distances 
from the centre of the disc r as a function of time: a) /i = 4 and 
a = 0.3; b) 12 = 1.8 and a = 0.3. 




time 

Figure 10. The change of the wave phase (pw(t) at fixed radius 
r = 0.67 (see figure 9a). There are two stages distinguished by 
their different angular velocities Qp = dcpw/dt: slow rotation with 
Qpi and a short phase of fast rotation with Qp2 > flpi 



A distinctive feature that we observe in our modelling 
is that the spiral pattern does not show a rigid-body rota- 
tion. Instead, at any radial distance in our models, its pos- 
sesses a two-pattern angular velocity, which switches from 
Qpi to Qp2 and back (see figure 10) twice during the orbital 
period. Most of the time the wave rotates more slowly (an- 
gular velocity r^pi), which corresponds to the case when 
the corotation radius rd locates at the periphery of the disc 
or beyond it. When the wave passes through the dark halo 
potential well, its angular rotation rises to Qp2 > ^pi, and 
the corotation radius rc2 moves inward down to (2 — 3) ■ rd- 

3.5 Models with the live halo 

To analyse the influence of possible interior substructure 
in the dark matter, e.g. of presence of s ubhaloes pre- 
dicte d by the standard cosmology theory (|Klvpin et al.l 
I1999I : iBullock et aD I2OOII : iRocha et all l2012h . we ran our 
A/'-body simulations of interaction between the triaxial live 
halo and the stellar collisionless disc. We performed the sim- 
ulations in several steps. First, we construct a quasi-steady 
triaxial live halo with a fixed stationary disc (a "frozen 
disc") inside. The halo is characterised by the parameter 
£~ 0.2 — 0.3 as described before, which is initially assumed 
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Figure 11. Comparison of the disc evolution in the models with 
live halo: the solid line is the case of a non-axisymmetric halo 
e = 0.1, the dashed line is a case of symmetric halo e = 0. Top 
panel: evolution of the Fourier harmonics m = 2. Bottom panel: 
evolution of the velocity dispersion. 
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Figure 12. The spiral pattern (left) and the halo isodenses in 
projection to the disc's plane (right) in our model with the live 
halo. 

to be constant along the radius. All disc particles were po- 
sitioned within the radius Ropt = 1.5, and the halo particles 
were extended to i?/^ = 5. The total number of particles was 
4 • 10^ for the halo and 2 • 10^ for the disc. The total halo 
mass in our modelling was Mh — 4. After several dynamical 
times the halo particles reached the state of relaxation and 
the initial numerical noise dissolved. Then we "unfreeze" the 
stellar disc particles instantly. Right after that, a two- arm 
spiral structure started forming, similar in its kinematics 
and morphological properties to that appeared in the case 
of rigid halo (see figure 12). 

It is important to consider the effects of the stellar disc 
feedback on the dark halo shape in projection to the disc 
plane for a wide range of the non-axisymmetry parameter e. 
In figure 12 we demonstrate the isodensity contours of the 
dark halo projected into the disc's plane at t = 12 (about 
4 rotation periods at Ropt)- We trace the evolution s of 
the halo at three regions separately: 1) the central region 
^1 ^ Td = 0.25, 2) the disc periphery at rs ^ 1, and the in- 
termediate region ri < r < rs (see figure 13). Up to time t = 



0.5 r — © — I — I — I 

o r<0.25 A r>l 




' — ^ — ^ — ^ — " — ^ — ^ — ^ — ^ — ^ — ^ — ^ — ^ — ^ — e — I 
4 8 12 

time 

Figure 13. Evolution of the e in the live halo at different radii. 
The line shows the best fitting for the central part of the disc. 

10, the amplitude of the spiral perturbations remains small 
and there is no noticeable stellar disc feedback on the halo 
shape. Some energy initially transfers from the triaxial halo 
particles to the disc, causing the spiral pattern generation. 
After that the central halo region becomes more axisym- 
metric, which means some energy is returned back to the 
halo due to the interaction with the spiral pattern. Local 
non-axisymmetric parameter in the intermediate zone also 
decreased down to ^ 0.15. The influence of spiral waves on 
the halo shape in the outer region remains insignificant. 

We executed several simulations with live and symmet- 
rical halo to investigate whether the disc and halo particle 
interactions can produce the disc heating. Such dynamical 
interaction could also generate some transient structures. 
Our modelling demonstrates the lack of the dynamical disc 
heating, at least in the marginally stable stellar disc in the 
presence of the live DM halo. In early moments of the mod- 
elling some oscillations of the disc velocity dispersion caused 
by relaxation of a small initial imbalance in particles distri- 
bution is observed (see the bottom panel in figure 11). After 
the relaxation is finished, the increasing of the velocity dis- 
persion stops. Transient structures forming in the disc have 
an insignificant amplitude < 1% in the case of the symmetric 
halo (see e.g. top panel in figure 11). Note that its ampli- 
tude is several times larger than in the model with rigid halo 
model. 



4 DISCUSSION 

Cosmological simulations show that dark matter plays sig- 
nificant role in formation and evolution of galaxies of all 
morphological types. In this paper we consider a possibil- 
ity of generation of a spiral pattern in the presence of non- 
axisymmetric dark halo. 

The key question related to the spiral structure in 
the o bserved galaxi es is whether the sp irals are tran- 
sient (|Sellwoodll201lh or long-lived features (|Bertin and LinI 
'1996'). In this paper we bring the existing evidences of non- 
spherical and, in general case, a triaxial dark matter distri- 
bution within and around the galaxies. Our numerical sim- 
ulations suggest that the long-lived spiral patterns in stellar 
discs may be generated not only by the internal mechanisms 
(as it is traditionally assumed) but also via gravitational 
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interaction between A triaxial dark halo and the galactic 
disc. Large-scale spirals emerge in our modelling even if the 
halo's triaxiality is as small as e ^ 0.02 — 0.05. We demon- 
strate that mechanism of the long-lived spiral pattern for- 
mation discussed above acts very effectively even in dynam- 
ically hot stellar discs (initial value of the Toomre param- 
eter Qt 2 — 4). Our preliminary simulations shows that 
the presence of a gase ous component does not change sig- 
nificantly this picture ([Khoperskov et aL I I2OI2I ). Note that 
similar mechanism can explain the existence of grand de- 
sign S2iral_Qattern dwarf early type galax- 
ies (|jerien et al.l I2OOOI : iLisker et al.l I2OO6I ) : we can assume 
that their spherical components are slightly triaxial. 

We also infer that the halo axes orientation can be found 
from kinematic studies of the spiral pattern. As we showed 
above, the pattern angular velocity cannot be described by a 
single constant value, as in the classical axisymmetric mod- 
els: at any radial distance it changes periodically twice for 
the period of rotation, as a result of nonaxisymmetry of the 
gravitational field of the halo (see figure 10). In principle, 
detailed measurements of the stellar velocity field can help 
detect the jumps of the phase velocity in the pattern (see 
figure 10), which spatially should point us at the major axis 
of the halo. The amplitude of the angular velocity pertur- 
bations can trace the deviation from the symmetry of the 
halo. The best candidates to this halo diagnostics are iso- 
lated non-barred spiral galaxies with small pitch angle of 
their spiral arms (^ 20°). Apparently the present day accu- 
racy of the methods of the Qp estimation makes tracing its 
variation along the radial and azimuthal direction a difficult 
problem. Note however that there are some observational 
evidences of multiplicity of the pattern angular velocity in 
the inner and in the ou ter regions of galactic discs (see f.e. 
iButa and Zhan j (|2009h V A combined estimation of the bar 
pattern sp eed in th e Milk y Way obtained from the gas dy- 
namics bv iGerhardI (|201lh gives — 52 zb 10 kms~^ kpc~^ 
and places the corotation radius near the end of the bar at 
3.5 — 5 kpc). On the other hand, the Milky Way spiral arms 
rotate much slower: observations constrain the speed of the 
spiral pattern in th e range of ^ 17 — 28kms~^kpc~^ 
(|Chakrabartvll2007f ). what corresponds to the corotation ra- 
dius position near the solar neighborhood 8 kpc). Thus, 
it is possible that the outer spiral structure in the Milky 
Way is relevant mostly to the triaxial halo-disc interaction, 
whereas the inner spiral arms are connected with the Milky 
Way bar. 

We demonstrate in figure 13 how the dark matter halo 
becomes more spherical in the centre due to the spiral struc- 
ture formation in the case of the initially triaxial halo. How- 
ever, we bring detailed consideration of effects of the an- 
gular momentum feedback from the disc to the halo be- 
yond the scope of this research. The dynamical interaction 
between the stellar disc and the DM halo occurs at reso- 
nan ces, and the result i s the increasing of the halo rota- 
tion lAthanassouTal (|2003l ). It is worth mentioning that the 
angular velocity of the bar in barred galaxies decreases in 
time due to dynamical frictio ns agains t the dark matter 
halo (see f.e. iDebattista and Sellwood (|2000l lV Numerical 
cosmological simulations reveal the possibility of variations 
of the axes ratios in the galactic haloes during their evo- 
lution. Here we note tha t para llel to the secular processes 
(e.g. see IDebattista etaP (|2006l )) the exchange by the angu- 



lar momentum and energy between the live halo and density 
waves it creates may also lead to the quasi-stationary peri- 
odical variation of parameters of the spiral arms (including 
their pitch angle) with time. The latter, in turn, may cause a 
variation of the n iorphological type o f a sp iral galaxy during 
its evolution fsee iMazzei and Curirl (|2003l )). 



5 CONCLUSIONS 

We conducted high resolution A/'-body simulations to study 
interactions between the stellar disc and triaxial massive 
dark halo in terms of the spiral structure formation in disc 
for different initial disc and halo parameters. Our results are 
most relevant to the galaxies with mass-geometric proper- 
ties similar to the Milky Way. The main conclusions are as 
follows: 

1. A triaxial dark matter halo can produce and maintain 
a long-lasting (at least a few Gyr) grand design spiral struc- 
ture in a stellar disc even if the disc is gra vitationally stable. 

2. A s in the case of gaseous discs (see iKhoperskov et all 
(|2012l )). the stellar discs respond to the dark halo asymmetry 
via the global spiral pattern formation, even if triaxiality in 
the halo is rather small (e — 0.02 — 0.05). 

3. The resulting amplitude of the spiral waves (m = 2) 
depends on the triaxiality s and is practically independent 
of the initial disc conditions. 

4. Morphology of the spiral structure in our models is 
rather similar to that formed in classical density wave the- 
ory. The growing spiral pattern does not lead to significant 
redistribution of mass in the disc in the radial and azimuthal 
directions in the case of a non- axisymmetric massive halo. 

5. Triaxial halo does not result to significant dynamic 
heating of the stellar disc (except, perhaps, the very initial 
stages of the interaction), in contrast to the case of a spiral 
pattern formed due to gravitational instability mechanisms. 

6. Dynamics of the spiral pattern cannot be described by 
a single pattern angular velocity. Instead, we observe ob- 
viously two pattern angular velocities at any fixed radius. 
Most of the time the spiral pattern at given radius rotates 
slowly with the angular velocity Qpi , but twice for it rotation 
cycle its angular velocity jumps up to Qp2 ^ (4 — 5)Qpi and 
then returns back. These jumps is caused by passing through 
the elongated potential well of the dark matter halo. 

7. In general, morphology and kinematics of the spiral 
pattern in a stellar disc are quite similar to those ob- 
tained in modelling of gaseous discs embedded in non- 
axisymmetrical dark matter haloes (|Bekki and FreemanI 
I2OOI iKhoperskov et al.ll2012l ) . 

8. There is a tendency for the halo to become more sym- 
metric in the disc plane in the models that take into ac- 
count self-consistent interaction between the disc and the 
dark halo. In this case a formation of spiral density waves 
remains effective in the outer regions of disc even when the 
inner regions of the halo have become axisymmetric. 
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